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We investigate the long-term dynamical evolution of star clusters in a steady tidal field produced by its 

parent galaxy. In this paper, we focus on the influence of mass profile of the parent galaxy. The previous 

\^ , studies were done with the simplification where the parent galaxy was expressed by point mass. We express 

different mass profiles of the parent galaxy by the tidal fields in which the ratios of the epicyclic frequency 

Q^l K to the angular velocity w are different. We compare the mass- loss timescale of star clusters whose tidal 

I ■ radii are identical but in parent galaxies with different mass profile, by means of orbits calculations in fixed 

O I cluster potential and N-body simulations. In this situation, a cluster rotates around the parent galaxy more 

■ ■ rapidly as the parent galaxy has shallower mass profile. We found that the mass-loss timescale increase 

^ ' 20% and 50% for the cases that the mass density profile of the parent galaxies are proportional to R~^ 

and R~^-^ where R is the distance from the galaxy center, compared to the point-mass case, in moderately 

strong tidal field. Counterintuitively, a cluster which rotates around the parent galaxy more rapidly has 

a longer lifetime. The increase of lifetime is due to the fact that the fraction occupied by regular-like 

orbit increases in shallower profile. Finally, we derive an evaluation formula for the mass-loss timescale of 

clusters. Our formula can explain a property of the population of the observed galactic globular clusters 

fT^ . that their half-mass radii become smaller as their distances from the galactic center become smaller. 
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Table 1. Relation among mass density profiles of galaxies, the angular velocities of the clusters, and the epicyclic frequency. 
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1. Introduction 

In a star cluster, stars exchange their energy with the other stars through their encounters, and the whole system 
approaches to thermodynamically relaxed state. This process is called two-body relaxation. In this process, some stars 
have high energy enough to escape from the cluster and the system loses its mass. Therefore, the mass-loss timescale 
of a cluster is thought to be proportional to its two-body relaxation time (Ambartsumian 1938; Spitzer 1940). 

However, this is not true for a cluster in an external tidal field, which was first shown by means of iV-body 
simulations in the Collaborative Experiment (Heggie et al. 1998). Fukushige and Heggie (2000) (hereafter, FH) 
showed that stars which have enough energy to escape from the cluster (hereafter, "potential escapers") can remain 
inside the cluster before finding exits on timescale comparable to the two-body relaxation time of the cluster, and that 
the timescale (hereafter, "escape time delay") complicates the scaling of mass-loss timescale with two-body relaxation 
time. Baumgardt (2001) and Baumgardt and Makino (2003) performed TV-body simulations in the external tidal field, 
and found the mass-loss timescale arc proportional to the half-mass relaxation time to the power of 3/4. Tanikawa and 
Fukushige (2005) (hereafter, paper I) also performed A^-body simulations, and found that the power depends on the 
strength of the tidal field and the power 3/4 shown by Baumgardt (2001) can be observed in the case of moderately 
strong tide. These findings mean that for treatment of mass loss of clusters in the external tidal field, it is desirable 
for orbits of stars to be fully resolved, and A^-body simulation is more appropriate than orbit-averaged method, such 
as Fokker-Planck method. 

This is a succeeding paper of Paper I. In this paper, wc investigate the effect of the mass profile of the parent galaxy 
which produces the tidal field. In most of previous simulations of star clusters with the external tidal field (Gicrsz, 
Heggie 1994; McMillan, Hut 1994; Fukushige, Heggie 1995; Baumgardt 2001; Paper I; Trenti ct al. 2007), the parent 
galaxy are simplified as point mass. 

If we consider the mass profile of the parent galaxy, the equations of motion of star clusters moving in the circular 
orbit around the galaxy are expressed, with tidal approximation, as 

-V$,.,:, (1) 
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where w is a angular velocity, k is a epicyclic frequency, ;y is a vertical frequency. The first term on the right-hand side 
in equation (1) is the Coriolis force, the second term is a combination of centrifugal and tidal forces, and the third 
term is gravitational acceleration from other particles in the cluster. For the point-mass galaxy, k = u. Also for the 
galaxy with flat rotation (p oc r^"^), n ~ 2a;. In Table 1, the approximated relations between mass profile and k/oj are 
summarized. The exact description will be given in section 2. 

In this paper, we compare dynamical evolution of the star cluster in the tidal fields of patent galaxies with four 
different k'^/uj'^{= 1,2,2.5,3), i.e. different mass profile. Throughout the paper, we set the initial tidal radii identically 
among clusters with different k/uj. The tidal radius and that at the initial time are, respectively, expressed as 

GM ^ ^/^ 



and 

rt,i = (GMi/(4c.2--«;2))i/3^ (3) 

where G is the gravitational constant, M and Afj are the cluster mass and that at the initial time. In the condition, 
the Coriolis force is relatively stronger for larger k} juP', i.e. shallower mass profile. Although several studies set 
parent galaxy to have mass density profile other than point mass (Portegies Zwart et al. 2001a; Portegies Zwart et 
al. 2002; Portegies Zwart et al. 2003; Moore 1996; Baumgardt and Makino 2003; Portegies Zwart et al. 2001b; Yim 
and Lee 2002; Dehnen et al. 2004; Capuzzo Dolcetta et al. 2005; Combes et al. 1999), they aim to study disruption of 
star clusters taking place in relatively short timescale. On contrary to them, wc investigate systematically the effects 
of the difference of mass density profile of parent galaxy on long-term mass-loss timescale of clusters . 
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The plan of this paper is as follows. In section 2, we give an exact description on the equation of motion in the tidal 
field. In section 3, we first perform orbital calculation in the fixed potential of star clusters and we compare the escape 
time delays of the potential escapcrs, which complicates the scaling with two-body relaxation timescale. In section 
4, we perform A'^-body simulation of clusters in the tidal field in order to study the infiucncc of mass density profiles 
of the parent galaxy on the mass- loss timescale of the star clusters. In section 5, we apply our simulation results for 
galactic globular clusters. Finally we summarize our paper in section 6. 

2. Equation of Motion 



In this section, we derive the equation of motion of stars in a cluster which moves circularly around a spherically 
symmetric parent galaxy with various mass density profiles. We set coordinate system which moves with the cluster. 
The origin of the coordinate system is set at the cluster center, the x axis is oriented in the opposite direction of the 
galactic center, and the z axis is along with the angular velocity vector of the cluster. Then, irrespective of the mass 
density profile of parent galaxy, we express the equation of motion of the cluster stars as 
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where r^ is the position vector of z-th star, J7 is the angular velocity vector of the cluster, Rg is the distance between 
the galactic center and the cluster center, and $G,i and $c,i are, respectively, the potential of the parent galaxy and 
the cluster at the position of the i-th star. The first and second terms in equation (4) are, respectively, the Coriolis 
and centrifugal forces on the i-th star. 

We consider the situation where the size of the cluster is sufficiently smaller than the distance from the galactic 
center to the cluster. Then, the galactic potential $g can be expanded in Taylor series around the cluster center, such 
as 
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where = (0,0,0) is the zero vector, the subscripts of indicate the galactic potential and its derivatives at 0. 
Since the galactic potential is spherically symmetric, {d^G/dy)o = (i9$g/9z)o = 0. Furthermore, {d^G/dx)o = Rgi^^, 
and {d^^Q/dy^)o = {d^^G/dz^)o = uj'^, where w is the angular speed of the cluster around the galactic center, i.e. 
rj = (0,0, w). Finally, equation (5) is rewritten as 
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Substituting equation (6) into equation (4), the equation of motion of the i-th star is expressed as 
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where 
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and K is so-called epicyclic frequency. 

The epicyclic frequency, k, and angular speed of the cluster, w, are determined by the mass density profile of the 
galaxy, po,-, and the distance between the galactic center and the cluster center, Rq,. The angular speed, w, is expressed 
as 
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where A^g(O) and pg(0) arc, respectively, the galactic mass and average mass density within the galactic radius, Rg. 
Since the mass density profile is spherically symmetric, the galactic potential, $Gi is expressed as 
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Differentiating partially equation (10) with respect to x, we obtain 



(10) 
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and therefore, 

dx' 
Substituting equation (12) into equation (8), we obtain the ratio of the epicycUc frequency to the angular speed, such 



==-2w^ + 47rGpG(0). (12) 
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where the second equality in equation (13) comes from equation (9). 

From equation (13), we can confirm the relation between the mass density profile of the galaxy and the ratio of k? 
to ixj"^ in table 1. Since pg{x) = in the case of pg{x) oc R~°°, we can see immediately k^ = w^. We consider the case: 

PG{x)^Cix + Rgy'{0<k<2). (14) 

Therefore, the average mass density within the galactic radius, Rg, is expressed as 

/° 4Tr{x + Rg)^pGdx 3 

where the second equality is obtained by using equation (14). Substituting equation (14) and (15 into equation (13), 
equation (13) is rewritten as 

— =A-klO<k<2). (16) 

The tidal radius is defined as the distance from the cluster center to a Lagrange point, and can be expressed by 
equation (2). In this paper, we investigate the cluster evolution by parameterizing the mass density profile of the 
galaxy, pQ, and the initial tidal radius given by equation (3), not the mass density profile of the galaxy, pG, and the 
distance between the galactic center and the cluster center, Rg. We regard the initial tidal radius, rt,i as one of the 
most important parameter for the mass-loss timescale of the cluster, since in paper I we found that the longer the 
mass-loss timescale of the cluster becomes, the larger the initial tidal radius is even if the mass density profile of the 
galaxy is equal. 

We consider the situation where the initial tidal radius, rt^i, i.e. (40;^ — k^), is kept constant. From equation (13), 
the angular speed, w, becomes larger, as the mass density profile of the galaxy becomes shallower, or pg(0)/pg(0) 
becomes smaller. In the equation of motion (7), the Coriolis force becomes larger as the mass density profile of the 
galaxy becomes shallower. 

As an example, we investigate the difference of the Coriolis forces between point-mass galaxy {pg{x) = 0,.t > — -Rg) 
and isothermal sphere galaxy {pg{x) = 3MG{0)/{4:TrRG{x + Rg)^),x > —Rg)- The values of K^juP' in point-mass 
galaxy and isothermal sphere galaxy are 1 and 2, respectively. Then, Coriolis force increases by 20 % from the case of 
point-mass galaxy to that of isothermal sphere galaxy. The Coriolis force changes angular momentum of cluster stars 
around cluster center. The difference of the Coriolis force may affect escape of cluster stars. 

3. Orbital Calculation of Potential Escaper in Fixed Potential 

In this section, we perform orbit calculations in the fixed and smooth cluster potential in order to investigate the 
escape time delay of the potential escapers. Due to their mutual encounters, stars gain enough energy to escape from 
the clusters and become the potential escaper. Such process takes place in two-body relaxation timescale. In the 
tidal field, the potential escapers can escape only through small apertures around the two Lagrangian points, and can 
remain inside the clusters before the escape on timescale comparable to two-body relaxation time (FH). Therefore, 
such escape time delay can influence on the escape timescale or the mass-loss timescale of the star clusters. 

3.1. Method 

We calculate the orbits of the potential escapers in a fixed and smoothed potential due to cluster stars plus the 
steady tidal field of the parent galaxy. The equation of motion of the potential escapers is described in equation (7). 
For the cluster potential, <I>c,i, we use Wq = 3 King model, where Wq indicates initial dimensionless central potential. 
In order to compute the acceleration from the King model, V^^i, in equation (7), we use cubic spline interpolation 
between grid data obtained by a numerical solution (e.g. Press et al. 1992). 
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Regardless of the mass density profile of the parent galaxy, we set the initial tidal radius of the cluster, rt^i to be 
equal to the radius beyond which the density is zero in the Wo = 3 King model, denoted by rkg. The value of rkg is 
3.13, where we use standard units (Heggie, Mathieu 1986), such that Mi = G = — 4£'c = 1, where Ec is the initial total 
energy within the cluster. When rtj = r^g = 3.13, the values of k and lo are determined for each given ratio k^ /uj'^, 
using equation (13) and (2). We performed the orbital calculations for 4 models in which the parent galaxy has mass 
density profiles of k'^/co^ = 1,2,2.5,3. The values of k^ and w^ are showed in table 2. The second right column shows 
representative mass density profiles of the parent galaxy in each ratio of k^/w^, such as single power law models. 

We determine the escape time delay as the duration between i = and the time when escape condition, \xi\ > rt, is 
satisfied. Note that the tidal radius, rt, remains the initial tidal radius, rt^j, during the calculation, since clusters do 
not lose the mass in this calculations. 

We calculate sets of orbits of the potential cscapcrs for several values of relative energy excess, i?pot,i = {Epot,i — 
-E'crit)/|£'crit|j whcrc -Bpot.i IS the energy of i-th potential cscaper, given by 

EpotAt) = 2 (^'x*^ + wyj^ + V.J) + $c,» ^2^^^"^- ^oj'^)x^ + 2^'^Zi^^ (17) 

where Wxj, "Uyi, and Wzi are x, y, and z components of the velocity of i-th particle, respectively. 

We integrate the orbits of the stars by means of adaptive stepsize control for a forth-ordcr Rungc-Kutta method 
(e.g. Press et al. 1992). The integration error in energy is about 9 x 10^^ at maximum {Epot.i = 0.03 at f = 3 x 10^ in 
standard units in the shallowest mass density profile of parent galaxy, k^ juP' ~ 3). 

3.2. Results 

In this subsection, we present the results of the orbital calculations of the potential escapers. In section 3.2.1, we 
investigate the difference of orbit families in planar orbits {z — plane) among parent galaxies with different mass 
density profiles. As the profile of the parent galaxy is shallower, more potential escapers move in regular-like orbit, 
and have larger escape time delay. In section 3.2.2, we show that this is also applicable to non-planar orbit families. 
In section 3.2.3, we investigate the distribution of the escape time delay in order to see whether the difference of orbit 
families affects the escape rate of potential escapers. 

3.2.1. Orbital families in planar orbits 

At first; our investigation is confined to the orbits of potential escapers on 2; = plane for simplicity of the analyses. 
If we give the potential escapers initial conditions where Zi = and Zi = 0, their orbits arc confined to z = plane 
forever shown in equation (7). 

In figure 1, we illustrate orbits and surfaces of section for 10 potential escapers with various escape time, te, for 
each model of k^/w^ = 1,2,2,5,3. We plot the orbit until the escape or t = 500 for non-escaper, and the surfaces of 
section until the escape or t = 10^ for non-escaper. Here and hereafter, potential escapers which are still in clusters 
when orbit calculations are terminated are called "non-escaper". We define the surface of section as {x,i) at which 
the potential escapers cross x axis in the direction of y> 0. 

At initial, we distribute 1000 potential escapers of i?pot,i = 0.03 uniformly in (x,y,x,y) phase space, and calculate 
the escape time delays. We sort the 1000 potential escapers in small order of the escape time delay. We chose fc-th 
potential escapers, where A: = 50+100(i — l)(i = 1,2,. . .,10) (hereafter, we call k{i)-th potential escapers "i-th potential 
escapers" ) , and show the orbit and surface of section of i-th potential escaper in i-th panel of each model in figure 1 
from left to right. As for the non-escapers, the selection are done at random. 

In the models of two steeper mass density profiles (k^/w^ = 1 and 2), the orbits become regular-like gradually 
from left to right. In both models, the 7 — 10th potential escapers looks like regular. However, in the models of two 
shallower mass density profile galaxies {k^ /oj^ = 2.5 and 3), the orbits are seen to become regular-like suddenly at the 
5th potential escapers. 

In all models, there are potential escapers whose surfaces of section are confined to small phase space, and they 
correspond to non-escapers. We call these orbits "regular orbits" . In the model of the shallow mass density profile 
(f? jitP' = 2.5), some orbits seem regular, but their surfaces of section occupy larger space than regular orbits, and 
eventually escape from the clusters, such as 6th and 7th potential escapers. We call these orbits "semi-regular orbits" . 

Wc consider the difference among the regular, semi-regular, and the other orbits. In figure 2 we put together surface 
of section of the 1st, 2nd, ..., 10th potential escapers in each model. The colors and shapes in surface of section 
are divided according to non- regular (white triangles) , semi-regular (white or gray squares) , and regular orbits (black 
circles), where we define semi- regular orbit as the orbit of potential escapers which have escape time delay, io > 100 
and non-regular orbit as escape time delay, io < 100. The surfaces of section of semi-regular orbits, whose escape time 
delay is ie > 100, are divided into two colors, white and gray. The difference of the colors shows not the difference 
among the potential escapers, but that of the time when the potential escapers cross x axis in the direction of y > 0. 
If the potential escapers cross x axis in the direction of y > less and more than 100 time units before they escape 
from the cluster, the colors of the surfaces of section are, respectively, white and gray. 
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The region in the surface of section occupied by the non-regular and semi-regular orbits is different from that of 
regular orbits in each model. The region in the surface of section for non-regular orbits is similar to that of semi-regular 
orbits just before escape, i.e. that of gray squares, and different from that of scmi-rcgular orbits not just before escape, 
i.e. that of white squares. From the above, we guess that the escape time delay of any potential escaper depends on 
their initial phase space. 

Figure 3 shows the escape time delay of potential escapers against the initial phase space {x,x). The uppermost 
4 panels show the potential escapers whose initial phase is y = 0, z = 0, y > 0, and i = 0. The points are colored 
according to the escape time delay, to, of the potential escapers. Black regions show non-escapers, and the lightest 
regions show the potential escapers of t^ < 10^. The escape time delay depends on the initial phase of the potential 
escapers. As the mass density profile of the parent galaxy is shallower from k^/cj^ = 1,2 to k} juj^ = 2.5, the phase 
volume for 10^ < io < 10^ becomes smaller, and that of t^ > 10^ becomes larger. Furthermore, as the mass density 
profile of the parent galaxy becomes shallower from k^/w^ = 2.5 to 3, the phase space for 10* <t^< 10^ is replaced 
by that of non-escapers. The phase space of tc < 10^ is similar in all models. 

In summary, we find that there is a larger fraction of regular-like orbital families in the parent galaxy with shallower 
mass density profiles. The difference in the orbital families among galaxy models brings the structure in initial phase 
space about their escape time delay shown in uppermost panel of figure 3. 

3.2.2. Orbital families in non-planar orbits 

Figure 3 shows the escape time delay against the initial phase {x,x) also for non-planar orbits. The points are 
colored according to the escape time delay of potential escapers. In each panels the escape time delay of the potential 
escapers whose initial positions and velocities of z-component arc in limited ranges are only plotted, and whose initial 
positions and velocities of y-component are y = 0, and y > 0. The z components increase from top panels to bottom 
panels. The structure in initial phase against the escape time delay seen in z=0 planer orbit can be also seen in other 
non-planer orbit, though the structure becomes less fine as the z components increases. 

3.2.3. Distribution of escape time delay 

Figure 4 shows fractions, i^pot, of potential escapers which do not escape at a given time t. The curves are drawn, 
from left to right, in order of large initial relative energy excess, Epot,i = 0.24,0.16,0.12,0.08,0.06,0.04 and 0.03. We 
spatially distribute them in the same way as primordial population in Wq — 3 King model immersed in a tidal field of 
parent galaxy. The vertical lines at the right side show the times when the integrations are terminated. A common 
feature in all models is that the potential escapers of larger relative energy excess, £'pot,i, have smaller escape time 
delay and smaller fraction of non-escapers. The feature in the case of the point mass (k^/w^ = 1) agrees with the 
result of FH. 

Figure 5 shows the fractions, Fpot, of potential escapers which do not escape at a given time t. We plot together 
the fractions, Fpot, of all models of k^/w^ in each relative energy excess, Epot,i = 0.03,0.08 and 0.24. The potential 
escapers in parent galaxy with shallower mass density profile escape more slowly, which agrees with the difference of 
orbit families described in section 3.2.1. 

We try to explain the difference of escape time delay among different i?pot,i and k^/w^ based on phase space flux of 
orbits of potential escapers. We define phase space flux near Lagrange points as amount of phase space which flows 
near Lagrange points in the outward direction with respect to the cluster per time unit. If potential escapers are 
equally present in every phase space, the number of escapers from the cluster should be proportional to the above 
phase space flux. In FH, which is the same as model of k^ juP' = 1, the phase space flux may explain the dependence 
of escape time delay on relative energy excess. 

The phase space flux, T , per unit energy is expressed (MacKay 1990) by 



J- = 2 / (5 



(j)' + -ix'^ + y'^ + z'^)-E 



x'dx'dy'dz'dy'dz', (18) 



where x' ,y' ,z' ,x' ,y' , and z' are coordinate and velocity when the origin is located in a Lagrange point which have 
X = rt, and 0' is obtained by expanding the effective potential, </>, around the Lagrange point to second order. and 
(j)' is respectively expressed by 

(b^-^{K,^-4iO^)x^ + ^LU^Z^ + '^,, (19) 

where $£ = —GM/rt, and 

^' - E„,, ^ ^{^ - Au;')x'^ - ^{k' - Au')y'^ - \{n' - hu?)z'\ (20) 

Since there arc two Lagrange points, we must add coefficient of 2 to integral. In equation (18), 0' is evaluated at 
x' = 0. We integrate equation (18), and obtain 
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^^^B=£^^. (21) 

A familiar calculation shows that the phase space volume per unit energy is expressed by 



V = 47r / ^/2{E - (j))(fr, (22) 



where the integration is done over the inside cluster. This does not depend sensitively on E in the vicinity of _E == £'ciit, 
and so we evaluate it by a Monte Carlo technique. The phase space fluxes and volumes per energy in Ep^t = 0.03 and 
the timescale of phase space flux, ip^ = V/J-, are summarized in table 3. 

Figure 6 shows the fraction, i^csc, of escapers which do not escape at a given time scaled by the timescale of the 
phase space flux, i/tph. We exclude non-escapers from figure 4. As seen in figure 6, the i^esc curves are in good 
agreement in three steeper mass density profiles (k^/oj^ = 1,2,2.5). The result in the case of point-mass galaxy agrees 
with that of FH. We can see that escape time delay, U, is proportional to the timescale of the phase space flux, 
tph (X {E — i^crit)"^ oc E~^ in the case of steeper mass density profile of parent galaxy (k^ juP' < 2.5). However, in the 
shallowest mass density profile {k} juP' = 3), the agreement becomes worse. 

In Figure 7, the fraction, Fosc of Epot,i = 0.03 for all models are plotted together against time scaled by the timescale 
of the phase space flux, i/iph- Agreement in the curves is not good, and escape time delay scaled by the timescale of 
the phase space flux, ie/iph, are longer as the mass density profile of the parent galaxy becomes shallower. Figure 7 
shows that the timescale of the phase space flux may not well explain the dependence of the escape time delay the on 
mass density profile of the parent galaxy. 

4. AT-body Simulations 

In this section, we present results of iV-body simulations of star clusters, and investigate how the dependence of 
escape time delay on mass density profile of parent galaxy influence the mass-loss timescale of clusters. In section 4.1, 
we describe the simulation method. In section 4.2, we show the results, and flnd that mass-loss timescale becomes 
larger as mass density profile of parent galaxy is shallower. 

4.I. Simulation method 

Wc investigate mass loss of star clusters in an external tidal field by means of A^-body simulations. The equation of 
motion is described in equation (7). In our simulations we use a softened gravitational potential, and the third term 
of equation (7) is expressed as 

-V<l>„ = - E 7, ,2' ,,3/2 ' ^23) 

where rrij is mass of j'-th particle and £ is a softening parameter. We set the softening parameter, e, as 1/32, where 
we use the standard units, A/j ~G ^ — 4i?c = 1, as section 2. 

We use Wo = 3 King's models (King 1966) to generate initial distribution of star clusters. We perform eight sets of 
simulations of star clusters with different initial tidal radii, rt^, and different mass density profiles of parent galaxies, 
iP juP ^ as summarized in table 2. The number of particles used for runs are A^ = 2* in the range described in table 
2. All particles have the same mass, m = Mi/N . We perform five runs whose realization of particle distribution are 
different for each A^ when A^ < 8192, and one run when A^ > 16384, except for model rt^i/rkg = 2.2 and k} juP — 3. For 
this model, five runs when A^ < 4096, and one run when A^ > 8192. 

The simulation code is the same as that used in Paper I. We perform numerical integrations of equation (1) using 
a leap-frog integration scheme with shared and constant timestep. The stepsize, Ai, is set to be as 1/64 in models 
''t,i/''kg = 1-0 and tP juP = 1,2,2.5, and as 1/128 in the other models. We used the Barnes-Hut tree algorithm (Barnes, 
Hut 1986) on GRAPE-5 (Kawai et al. 2000), a special-purpose computer designed to accelerate A^-body simulations. 
We use only the dipole expansion and the opening parameter Q = 0.5. It spend about 500 CPU hours completing the 
longest run, N = 32768 in models which have larger tidal radii (rt,i/rkg = 2.2) and the steepest mass density profile 
of parent galaxy (jp juP = 1). For smaller A simulations, the force calculation is done by direct summations (when 
A £ 4096) and on host computer (without GRAPE-5, when A^ 5^ 512). 

Contrary to the standard star cluster simulations, our simulation uses a softened gravitational potential and a 
leap-flog integrator with a relatively large stepsize. Also, the force calculation is performed with the tree algorithm. 
We adopt these approaches, since we partly use mass-loss timescale obtained in paper I, which modeled mass density 
profile as the point-mass one {hP /uP = 1). As discussed in paper I, the approaches that we adopt do not influence the 
results concerning escape from the cluster. 
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4-2. Results 

Figure 8 shows evolutions of the total mass for all models. The curves indicate the decrease in mass of the cluster 
defined by a tidal boundary. We define geometrically the cluster member as all stars within the tidal radius from 
the center of mass of the cluster, which is expressed as equation (2). Since M depends on rt, itself, some iteration 
is usually required. We remove stars when they escape far enough (more than 4096 length units from the coordinate 
origin) . 

Figure 9 shows the mass-loss timescale of the clusters as a function of the initial half-mass relaxation time, irh.i in 
the upper panels. The lower panels show the mass-loss timescale of the clusters scaled by those in the parent galaxy 
with the steepest mass density profile (k^ juP' = 1). The mass-loss timescale is here defined as the time when 50% of 
the initial total mass is lost. Since we perform several runs for each model, we use the means of these runs as the 
mass- loss timescale, imioss- Table 4 gives the maximum deviations from the means among these runs. Since we use 
the potential softening and fix the softening parameters for all models, the initial half-mass relaxation time, irh,i: is 
expressed as 

Nrv. -3/^ 

irh,i = 0.138 — — -^^ , (24) 

Mii/2^i/2ln(0.079rha/£) 

where rh.i is the initial half-mass radius. The Coulomb logarithm. A, is expressed as the ratio Pmax/Pmin, where Pmax 
and Pniin a-re, respectively, the maximum and minimum impact parameters. Since the maximum impact parameter, 
Pmax, is about the system size, it is expressed as Pmax = o'^h.i- The minimum impact parameter, Pmin, is determined by 
the softening parameter, e. Therefore, Pmin = &£• A ratio of a/h~ 0.079 is taken so that the evolutions of the minimum 
potential in a cluster would have no difference between runs for three different softening parameters, e, as can be 
seen in the lower panel of figure 15 of paper I. Thus, the initial half- mass relaxation timescale, irh,i, was estimated as 
irh,i = 1200 X (A^/8192). The initial half-mass radii in the standard unit are rh,i = 0.84. 

Even if clusters have the same initial half-mass relaxation time, the clusters lose their mass more slowly in the 
parent galaxy with shallower profile. For models in a moderately strong tidal field (rt^i/r^g = 1.0), clusters in parent 
galaxy with shallower mass density profile (k^/w^ = 2,2.5,3) have larger mass-loss timescale than those in the parent 
galaxy with the steepest mass density profile (k^/w^ = 1) by 20 %, 50 %, and factor of 2.5 respectively. On the other 
hand, the largest difference of mass-loss timescale is at most 60 % for the case in weak tidal field (rt^/y'kg = 2.2). As 
the tidal field becomes stronger, sensitivity of mass-loss timescale to mass density profile of parent galaxy is enhanced. 

We consider possibilities other than escape time delay for the difference of mass-loss timescale among clusters in 
different parent galaxy. At first, we check whether the difference of mass-loss timescale is affected by tidal lock. Since 
clusters set so far rotate around itself at angular velocity (0,0, w) with respect to inertia frame, their rotation speeds 
around itself are different among clusters in parent galaxy with different mass density profiles. 

We perform iV-body simulations of clusters which rotate solidly at angular vector (0,0,— w) in the same reference 
frame described in section 2. The clusters have Wq = 3 King profile and are in a moderate tidal radius, rt,i/rkg = 1.0. 
These clusters do not rotate with respect to inertia frame. Figure 10 shows mass-loss timescale of clusters tidally 
unlocked, and mass-loss timescale scaled by mass-loss timescale of clusters tidally unlocked in the parent galaxy with 
the steepest mass density profile (k^/w^ = !)• 

Mass-loss timescale and the increase from steeper mass density profile of parent galaxy to shallower one are not 
distinguishable among clusters tidally locked and unlocked. Therefore, tidal lock is not responsible for the difference 
of mass-loss timescale among clusters in parent galaxy with different mass density profile. 

Next, we check whether two-body relaxation, which induces the escape of stars, is infiuenced by the difference of 
mass density profile of parent galaxy. However, we find that two-body relaxation docs not change. Figure 11 shows the 
change in energy of an individual star for clusters in parent galaxy with mass density profile k^/w^ = 1,2,2.5,3, where 
''t,i/''kg = 1-0, and the clusters have N = 1024. In figure 11 we plot {Smax.ijmod, the median value of the maximum 
energy records at a given time, ti: 

Smax,»(ii)-max{^,(f)}, (25) 

t<ti 

E^{t) = -(«x,' + t'y»' + vj) + $c,, - 2 ('*' " 4cj2)2;,2 + _z,2, (26) 

of 50 randomly selected particles. We did not update -Emax,i(''i) after stars escaped from the cluster. The representative 
value, {i?max,i}med, Hiay tracc the energy acquired by two-body relaxations, on the average. When stars obtain enough 
energy, they escape from cluster eventually. The dashed lines indicate the escape energy of each model, which is 
potential at Lagrange points and expressed by 

'iGM 

-E-crit — -z . (27) 

2rt 
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As shown in figure 11, steep increase of {-Emax/ilmcd from —0.8 to —0.6 are similar among all clusters. The increases 
of {i?max,i}med at i < 400 are very similar among parent galaxy with steeper mass density profile (k^/w^ = 1,2,2.5). 
The facts above indicates that the relaxation processes in all models occur on almost the same timescale. The deviation 
of the shallowest model [k^ /uP' = 3) from other models after t = 100 is due to the much slower mass loss than those 
of other clusters. 

5. Evaluation Formula and Implication for Observed Galactic Globular Clusters 

In this section, we formulate mass-loss timescale of globular clusters, using our simulation results. From the formula, 
we obtain properties of globular clusters which survive over the Hubble time. We compare the properties with those 
of the observed galactic globular clusters, and examine the validity of our formula for the mass-loss timescale. 

5.1. A formula of mass-loss tim,escale 

We derive an evaluation formula for the mass-loss timescale of a globular cluster using the result of Paper I and this 
paper that the mass loss timescale is proportional to the half-mass relaxation time, imioss ex ij^j^, and a correction factor 
for the difference of the mass-loss timescale among mass profiles of parent galaxies. We write the mass-loss timescale 
of a globular cluster, imiossi a-s 

tn.loss{n,^Vu;^)^t^.loss,t,,^ = to{rt) i^j [l + /(A^V^')3(^t)] -, (28) 

where 

n=rt.i/ru,i. (29) 

The product of the first and second factors in the right-hand side of equation (28) represents the mass-loss timescale 
of a cluster in a tidal field produced by a point-mass galaxy. The timescale imioss,trh=to(^t) is mass-loss timescale 
of a cluster whose initial half-mass relaxation time iih.i is a given time ig- The third factor in the right-hand side 
of equation (28) represents the increment from mass-loss timescale of clusters in a point-mass galaxy to those with 
shallower galaxy. 

We construct fitting formulae for four functions in equation (28), tmioss,t,h=to(^t), a;(ft), gi^t), and fi/PjuP'), using 
our simulation results summarized in table 5. For the function tmioss.t,.h=to(^t)! we use those of clusters whose initial 
half-mass relaxation time is 1000 A^-body time unit, iiooo(^t) = ^mioss,trh=iooot„u(^t)/inu, where tnu is the iV-body time 
unit. The simulation results of the mass-loss timescale tiooo(^t) are shown in the third column of table 5. The fitting 
formula of iiooo(^t) is given as 

iiooo(n) « tiooo,fit(^t) = 2 X \{f-^-±—^. (30) 

As seen in the upper panel of figure 12, equation (30) is well fit to the mass-loss timescale, iiooo(^t), over an intermediate 
range of ft = 3 - 20. Since a cluster immediately loses its mass in the strong limit of a tidal field, iiooo(u) — >■ in the 
limit of ft — 5- 0. Since mass- loss timescale of a cluster does not depend on the strength of the tidal field in the weak 
limit, iiooo(^t) should asymptotically approach to a constant value. 

For the function a;(ft), we use those at the initial half-mass relaxation time tj-h.i is lOOOtnu, which is shown in the 
fom'th column of table 5. The fitting formula of the function a; (ft) is given as 

.(ft)«.tit(ft) = ^4^^^. (31) 

As seen in the lower panel of figure 12, equation (31) is well fit to the logarithmic slope, a;(ft), over an intermediate 
range of ?~t. Since mass-loss timescale of a cluster does not depend on its initial half-mass relaxation time in the strong 
limit of a tidal field, a; (ft ) ^- in the limit of ft — ?> 0. Since mass-loss timescale of a cluster should be proportional to 
its initial half-mass relaxation time in the weak limit of a tidal field, a;(ft) asymptotically approaches to unity. 

We construct of the fitting formula for the increment from mass-loss timescale of a cluster in a point-mass galaxy 
to that in a shallower galaxy, f{K^/ijj'^)g{ft), as follows. We simplify the dependence of the increment on ft, such that 
the increment /(K^/u;^)f/(ft) in a strong tidal field is the same as that for ft = 3.9, and the increment f{K^/u)'^)g{ft) in 
a weak tidal field is zero. These are shown in the fifth column of table 5. The fitting formula of the function /(k^/oj^) 
is given as 

/ (.V-^) « /.t (.V-^) = ^, /' ,,. , (32) 

16 — (k^/w^) 

and we adopt a step function for the fitting function of g{ft) as 
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As seen in figure 13, equation (32) is well fit to the increment obtained by our simulation of k^/oj'^ = 2,2.5,3. The 
increment, 1 + /{k'^/uj'^), approaches to infinity at k^ juj^ —5- 4. 

5.2. Survivability of galactic globular clusters 

Using equation (28), we obtain parameters of a cluster which survives over the Hubble time, 10^" years. We compare 
the parameters with those of the galactic globular clusters, and examine the validity of our formula for the mass-loss 
timescale. For obtaining the mass- loss timescale of a cluster from equation (28), the following set of parameters is 
required: the initial tidal radius, rt,i, the initial half- mass radius, rh^j, the initial half- mass relaxation time, irh,i, and 
the mass profile of the parent galaxy, k^/uj'^. 

Here we rewrite equation (28) using observable parameters of galactic globular clusters, the initial total mass of the 
cluster. Mi, the average mass of the cluster stars, to, the initial half-mass radius, rh,i, the mass profile of its parent 
galaxy, k^/w^, the circular velocity of the cluster, Vc, and the distance of the cluster from the center of the parent 
galaxy, Rq. Adopting tg = lOOOfnu, equation (28) is rewritten as 



tmlossift,!^'^/^^) = tiooQ{ft)t„ 



I iih,i 



V iooo<„ 



x{ft) 



[l + f[K^lu:-^)g[r,)] 



(34) 



Substituting the four fitting formulae (30), (31), (32), and (33) into iiooo('^t), a;(ft), /(k^/oj^), and g{ft) in equation 
(34) respectively, we obtain an evaluation formula for the mass-loss timescale of clusters as 



[iloss 



{h,K'/uj'')KlxlQ\Y] 



(■1000 



,fit(n) 



iiooo,fit(rt =4) 



A/i 



2 X lO-^M 



G.3A/, 



ni,i 
3pc 

a;tit(n) 



3/2 



Ms 



2 X lO-^Afm 



-1/2 



[l + /fit(A^Vw')5fit(rt)], 



where ft is given by 



h = 3.5 



4-^7^2 



-1/3 



M; 



2 X 104 Aff. 



1/3 



/ R, 



\lkpc 



2/3 



200km/s 



-2/3 



Hi 
3pc 



(35) 



(36) 



Here we use the iV-body time unit is tnu ~ f\!i jG^I'^M-^ ' , and the initial tidal radius and half-mass relaxation time 
are, respectively, expressed as equation (3) and (24). 

We set our Galaxy as a model of the parent galaxies, and set to have the constant circular velocity of 200 km/s, 
such that k} jixP' = 2 and Vc = 200 km/s. The average mass of the cluster stars is set to be O.3A'/0. 

Solid curves in figure 14 shows the mass- loss timescale of clusters as a function of their initial half- mass radii, rh,i. 
The initial total mass of the clusters, Afj, are indicated in each panel. The distances of the clusters from the center 
of the parent galaxy, i?G, are described beside each curve. The dotted lines show twenty times the initial half-mass 
relaxation time. The dashed lines indicate the Hubble time, 10^° years. The kink of each curve around the peak of 
the mass- loss timescale is due to the discontinuity of the step function, g{ft). 

As seen in figure 14, the mass-loss timescale of clusters is proportional to their initial half-mass relaxation time, if 
their half-mass radii are sufficiently small. When the half-mass radii are larger than a critical half-mass radius, the 
mass-loss timescale of clusters becomes smaller as their half-mass radii become larger. This is because the tidal field 
produced by their parent galaxy strongly affects the mass loss of the clusters. The critical half-mass radius is larger 
as the distance of the cluster from the center of the parent galaxy is larger, i.e. the tidal field is weaker. 

Figure 15 shows the properties of the galactic globular clusters: their half-mass radii and their distances from the 
galactic center, which arc obtained from Harris (1996). We divide the galactic globular clusters by their luminosities: 
absolute visual magnitude M^ < — 7 (filled circles) and M^ > — 7 (open circles). The mass of a cluster with absolute 
visual magnitude Afv = —7 is about IQt'Mq, since the solar absolute visual magnitude is Afv ~ 5, and the mass-to- 
light ratio of the globular clusters is about 2. The dashed curves show parameters of clusters whose lifetimes are the 
Hubble time, derived from our evaluation formula (35). The values beside the dashed curves indicate the masses of 
the clusters. 

Both clusters with absolute visual magnitude My < — 7 and A/v > — 7 are absent above the upper limit of the initial 
half-mass radius in the case of a cluster mass 10^ Af©. The clusters with absolute visual magnitude A/v < ~7 may be 
absent from the beginning. We can say that our upper limit of the initial half-mass radii is in a good agreement with 
those of the galactic globular clusters in a wide range of the distance from the galactic center, Rq = 1 - 100 kpc. 

The globular clusters with half-mass radii less than their lower limit in the case of cluster mass lO^Af© (1 pc) are 
absent. This is consistent with the fact that most globular clusters have mass less than IO^ttiq. However, many 
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clusters with absolute visual magnitude M^ > —7, whose masses are less than lO^A/0, have half-mass radii less than 
their lower limit in the case of cluster mass lO^M©, ~ 2 pc. We may interpret that these clusters initially have larger 
mass than 10^A/q, suffer mass loss strongly, and become small, as seen nowadays. 

6. Summary 

We investigate mass loss of star clusters in a tidal field of parent galaxy with different mass density profiles. We find 
that orbit families are different among parent galaxy with different mass density profiles. A fraction of regular orbits 
increases in the orbit family as the mass density profile of the parent galaxy is shallower. The tendency is strong in 
nearly planar orbits. The difference of the orbit families affects the distribution of escape time delay in each parent 
galaxy model. The escape time delay becomes larger as the mass density profile of the parent galaxy is shallower. 

We conform that the difference of escape time delay affects mass-loss timescale of clusters with two kinds of size of 
tidal radii. The mass-loss timescale of clusters becomes larger as the mass density profile of the parent galaxy becomes 
shallower. In a moderate tidal radius, ^t^i/rkg = 1.0, the clusters in parent galaxy with the three shallower mass 
density profiles (k^/o;^ = 2,2.5,3) have larger mass-loss timescale than those in the parent galaxy with the steepest 
one (k^/cj^ = 1) by 20 %, 50 % and a factor of 2.5, respectively. In the case of larger tidal radius, rt,i/rkg = 2.2, the 
difference is smaller and the largest is 60 %, which is comparison between N = 128 clusters in parent galaxy with the 
steepest and the shallowest mass density profiles. 

The important result is that orbital family is attributed to the difference of the mass-loss timcscales among star 
clusters in parent galaxy with different mass density profiles. Even if potential escapers occupy the same phase 
space among different parent galaxy models, the orbits and the escape time delays may be different by an order of 
magnitudes. Thus, the difference of mass-loss timescale among different parent galaxies is not represented by orbit- 
averaged method, such as Fokker-Planck simulations. Orbit-averaged methods use escape criterion based on phase 
space of cluster stars, for example energy and angular momentum. 

We derive formula for the mass- loss timescale of clusters, as seen in equation (28). Our formula can explain a 
property of the population of the galactic globular cluster that their half-mass radii become smaller as their distances 
from the galactic center become smaller. 
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Table 2. Initial models. 



n.i 


rt,i/ri,g(Wo = 3) 


k'/u' 


K 


co'' 


p{R) 


N 


3.13 


1.0 


1 


0.011 


0.011 


p(xR-°° 


128- 


- 131072 


3.13 


1.0 


2 


0.033 


0.016 


p(xR-^ 


128- 


- 131072 


3.13 


1.0 


2.5 


0.054 


0.022 


p(xR~^-^ 


128 


- 32768 


3.13 


1.0 


3 


0.098 


0.033 


pocR-^ 


128 


- 32768 


6.97 


2.2 


1 


0.00098 


0.00098 


pozR-°° 


128- 


-131072 


6.97 


2.2 


2 


0.0029 


0.0015 


p(xR-^ 


128- 


-131072 


6.97 


2.2 


2.5 


0.0049 


0.0020 


p(xR-^-^ 


128 


- 32768 


6.97 


2.2 


3 


0.0088 


0.0029 


p(xR-^ 


128 


- 32768 



Table 3. Case of r't i/ri;g = 1.0. Phase space flux per unit energy, J^, phase space volume per unit energy, V in Epot,i = 0.03 and 
timescale of phase space flux, tpij. 



'2 I, :i 



K^/w 



T 



V ^ph 



1 0.108 317 2940 

2 0.102 311 3050 
2.5 0.0969 305 3150 

3 0.0885 296 3340 



Table 4. Maximum deviation of tmioss for ^ = 8192. 













(n.i/'^kg 


,k'/l,'') 










(1.0,1) 


(1.0,2) 


(1.0,2.5) 


(1.0,3) 


(2.2,1) 


(2.2,2) 


(2.2,2.5) 


(2.2,3) 
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32 


32 


39 


42 


40 


110 


71 


154 




256 


20 


22 


22 


39 


65 


115 


104 


86 




512 


19 


13 


100 


63 


90 


229 


61 


115 


N 


1024 


108 


67 


42 


105 


38 


189 


101 


182 




2048 


40 


54 


39 


42 


230 


130 


99 


158 




4096 


81 


53 


58 


238 


194 


248 


193 


227 




8192 


166 


58 


130 


133 


118 


362 


122 


- 



Table 5. Our simulation results used for our fltting formula 
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Paper I 
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Paper I 
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Paper I 
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Fig. 2. For each model, surfaces of section in figure 1 are put together. Surface of section is divided into four kinds, regular orbit 
(black circle), non-regular orbit (triangle), and semi-regular orbit (gray and white squares). The colors of the squares are divided, 
whether the potential cscapers cross x axis in the direction of y > loss (white) or more (gray) than 100 time units before escape. 
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Fig. 3. Escape time delay of potential eseapers as a function of their initial phase space {x,x). The color of initial phase space in 
each panel becomes deep from gray to black in order of escape time delay, te < 10^, 10^ < to < 10"^, 10^ < to < 10'', 10* < te < 10®, 
and non-escaper. The conditions about \z\ and \z\ show initial phase space of potential eseapers, and their (2,2) is distributed in 
the range at random. 
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Fig. 4. Distribution of escape time delay of potential escapers. Mass density profile of parent galaxy becomes shal- 
lower from left panel to right one. In each panel, the potential escapers have, from left to right, relative energy excess 
Epot.i = 0.24,0.16,0.12,0.08,0.06,0.04 and 0.03. 
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Fig. 5. Distribution of escape time delay of potential escapers. Relative energy excesses of potential escapers are 0.03,0.08 and 
0.24, from left panel to right one. 
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Fig. 6. Fractions of escapers which do not escape at a given time scaled by timescale of phase space flux, ipj,. Mass density profile 
of parent galaxy becomes shallower from left panel to right one. 
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Fig. 7. Fraction of cscapcrs which do not escape at a given time scaled by timescale of phase space flux, tp^. The cscapers have 
the relative energy excess Spot,; = 0.03. 
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Fig. 8. Evolution of cluster masses. We show an initial tidal radius of each cluster, rt/rj^g, and steepness of mass density profile of 
each galaxy, k? /ui^ , outside each panel. The numbers under the curves indicate the numbers of particles, A'^. 



18 



A.Tanikawa and T.Fukushige PASJ 



[Vol. , 



^t,/^kg=l-0 



^tAkg=2.2 




2.5 - 



II 



^ 2h 



— 1 1 — I — I I I I I I 1 1 — I — I I I I 1 1 1 1 — I — I I I I I I r- 



K-Vco2=3 




1 B .a..a.i □. a.a..D h. □..d q - 



— 1 — I — I — I — I I 1 1 1 



— 1 — I — I — I — I I 1 1 1 1 — I — I — I — I I I 

2.5 ^ - 
2 ■ 

1 ^ ■ 



10^ 



10' 



10' 




?rh,i 



^rh,i 



Fig. 9. Upper panels: mass-loss timescale of clusters as a function of initial half-mass relaxation time. Lower panel; mass-loss 
timescale of clusters scaled by mass-loss timescale of clusters in the parent galaxy with the steepest mass density profile {k?" jiJ^ = 1) 
as a function of the initial half-mass relaxation timescale. 
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Fig. 10. Mass-loss timescale of clusters tidally unlocked as a function of initial half-mass relaxation time (left panel). Mass-loss 
timescale of clusters tidally unlocked scaled by mass-loss timescale in which parent galaxy have the steepest mass density profile. 
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Fig. 11. Time evolution of the energy of an individual star (the definition is in the text) for clusters (A' = 1024) in a moderate 
tidal radius, rt,i/rkg = 1.0. The dotted lines indicates the time evolution of the escape energy of clusters in parent galaxy writh mass 
density profile in order of steepness from top to bottom. 
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Fig. 12. Upper panel: mass-loss timescale of a cluster whose initial half-mass relaxation time is 1000 A'^-body unit as a function 
of the ratio of the initial tidal radius to the half-mass radius, ft. Lower panel: logarithmic slope, x [= rfln(injioi3s)/dln(irij ;)], as a 
function of ft. The filled circles indicate the data points of our simulation results in paper I, and the solid curves show the fitting 
formula. 
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Fig. 13. Increment of mass-loss timescale, 1 + f{K^ /uj^), as a function of mass density profile (k^/ui^) with respect to the steepest 
one (k^/oj^ = 1). The filled circles indicate the data points of our simulation results in the lower left panel of figure 9, and the solid 
curves show the fitting formula. 
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Fig. 14. Mass-loss timescale of clusters as a function of their half-mass radii, rh,i. The initial cluster mass, Mi, is written on the 
top left side of each panel. The value beside each curve indicates the distance of the cluster from the galaxy center, Rg. The dotted 
lines show twenty times their initial half-mass relaxation time, i,.h,i. The dashed lines indicate the Hubble time, lO'^'' years. 
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Fig. 15. Half-mass radii of the galactic globular clusters, and their distance from the galactic center, obtained from Harris (1996). 
The clusters are divided by their luminosities: absolute visual magnitude A/v < —7 (filled circles) and Mv > —7 (open circles). 



